{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# MBQC 模型下求解多项式组合优化问题\n",
    "\n",
    "<em> Copyright (c) 2021 Institute for Quantum Computing, Baidu Inc. All Rights Reserved. </em>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "在[基于测量的量子近似优化算法](QAOA_CN.ipynb)中，我们介绍了多项式无约束布尔优化问题（polynomial unconstrained boolean optimization problem, PUBO），并提出了**基于测量的量子近似优化算法（MB-QAOA）** 对该类问题的求解方式，感兴趣的读者可以参阅前面的内容。这里我们给出两个示例作为 MB-QAOA 的实战演示：一个具体的 PUBO 问题和一个最大割问题。\n",
    "\n",
    "## 示例：PUBO 问题\n",
    "\n",
    "我们首先简单回顾一下什么是 PUBO 问题。给定一个变量为 $x = \\{x_1,\\cdots,x_n\\}$ 的 $n$ 元多项式，\n",
    "\n",
    "$$\n",
    "C(x) = \\sum_{\\lambda \\in \\Lambda } \\alpha_{\\lambda} \\prod_{j \\in \\lambda} x_j,\\tag{1}\n",
    "$$\n",
    "\n",
    "其中每个变量 $x_i \\in \\{0,1\\}$，$\\prod_{j \\in \\lambda} x_j$ 为一个单项式，$\\lambda \\subseteq [n]:= \\{1, 2, ..., n\\}$ 为一个指标集，$\\Lambda$ 为指标集的集合，$\\alpha_\\lambda$ 为每个单项式对应的实系数。在 PUBO 问题中，$C(x)$ 称为目标多项式。PUBO 问题就是寻找一组最优解 $x^* = \\{x_1^*, x_2^*, ..., x_n^*\\} $ 使得目标多项式的取值最大，即\n",
    "\n",
    "$$\n",
    "x^* =\\underset{x}{\\text{argmax}} \\ C(x).\\tag{2}\n",
    "$$\n",
    "\n",
    "代码实现上，我们规定输入多项式的标准形式为一个列表，列表中的第一个元素为多项式的最大变量数，列表的第二个元素为一个字典，字典的键为每个单项式涉及的变量（常数则为 ‘cons’），变量之间用逗号隔开，值为该单项式的系数。比如，我们希望输入一个三元多项式 $x_1 + x_2 - x_3 + x_1 x_2 - x_1 x_2 x_3 + 0.5$，应输入如下代码："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 变量数\n",
    "var_num = 3\n",
    "# 以字典格式存储多项式的各个单项式\n",
    "poly_dict = {'x_1': 1, 'x_2': 1, 'x_3': -1, 'x_1,x_2': 1, 'x_1,x_2,x_3': -1, 'cons':0.5}\n",
    "# 构造符合约定的多项式\n",
    "polynomial = [var_num, poly_dict]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**注意：** 由于我们考虑的变量为布尔值，所以每个单项式中变量指数最大为 1，即输入的多项式字典的键中，每个变量最多出现一次。例如：{'x_1,x_1,x_2': 1} 为不符合规范的输入。我们同样约定自变量的下角标从 \"1\" 开始，与数学习惯相一致。并且为了统一输入规范，自变量需要连续编号，不能出现输入多项式为 $x_1x_2 + x_6$ 的情况，应该改为 $x_1x_2 + x_3$，否则我们会以报错的形式提醒输入不符合规范。\n",
    "\n",
    "代码实现上，我们在 `pubo` 中定义了 `is_poly_valid` 函数，用于检查输入多项式的合法性。如果合法，将会打印 \"The polynomial is valid.\" 字样，否则会报错。\n",
    "\n",
    "```python\n",
    "from paddle_quantum.mbqc.QAOA.pubo import is_poly_valid\n",
    "```\n",
    "我们同样定义了一个函数 `random_poly` 根据变量个数随机生成一个布尔型多项式。\n",
    "```python\n",
    "from paddle_quantum.mbqc.QAOA.pubo import random_poly\n",
    "```\n",
    "**注意：** 随机生成的多项式并不一定是合法的！因此我们需要在计算之前检查多项式的合法性。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 求解 PUBO 问题的代码实现"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 引入记时模块\n",
    "from time import perf_counter\n",
    "# 引入符号计算模块\n",
    "from sympy import symbols\n",
    "#引入 paddle 相关模块\n",
    "from paddle import seed, optimizer\n",
    "# 引入 pubo 模块\n",
    "from paddle_quantum.mbqc.QAOA.pubo import dict_to_symbol, is_poly_valid, brute_force_search\n",
    "# 引入 qaoa 模块 \n",
    "from paddle_quantum.mbqc.QAOA.qaoa import MBQC_QAOA_Net, get_solution_string"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "利用 MB-QAOA 算法，我们可以定义 PUBO 主函数 ``mbqc_pubo``，输入目标函数和 QAOA 算法深度，最终返回最优解对应的比特串。在主函数中，**最核心的是 MBQC_QAOA_Net 类**，其集成了 MB-QAOA 以及优化网络，感兴趣的读者可以返回到 [基于测量的量子近似优化算法](QAOA_CN.ipynb) 中查看，这里我们直接调用。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 定义 MBQC 模型求解 PUBO 问题的主函数\n",
    "def mbqc_pubo(OBJ_POLY, DEPTH, SEED, LR, ITR, EPOCH, SHOTS=1024):\n",
    "    \n",
    "    # 将字典类型的目标多项式转化为 sympy 能处理的符号多项式\n",
    "    obj_poly = dict_to_symbol(OBJ_POLY)\n",
    "    var_num, poly_symbol = obj_poly\n",
    "\n",
    "    # 打印当前 QAOA 算法电路深度\n",
    "    print(\"QAOA 算法深度为：\", DEPTH)\n",
    "\n",
    "    # 开始计时器\n",
    "    start_time = perf_counter()\n",
    "    \n",
    "    # 初始化一个 MB-QAOA 训练网络\n",
    "    seed(SEED)\n",
    "    mbqc_net = MBQC_QAOA_Net(DEPTH)\n",
    "    \n",
    "    # 选择 Adams 优化器\n",
    "    opt = optimizer.Adam(learning_rate=LR, parameters=mbqc_net.parameters())\n",
    "\n",
    "    # 开始训练\n",
    "    for epoch in range(EPOCH):\n",
    "        # 每次迭代后更新训练参数\n",
    "        for itr in range(1, ITR + 1):\n",
    "            # 用网络进行训练并返回损失函数值和当前输出量子态\n",
    "            loss, vector_out = mbqc_net(poly=obj_poly)\n",
    "            # 根据损失函数值优化参数\n",
    "            loss.backward()\n",
    "            opt.minimize(loss)\n",
    "            opt.clear_grad()\n",
    "            if itr % 10 == 0:\n",
    "                print(\"iter:\", itr, \"  loss_MBQC:\", \"%.4f\" % loss.numpy())\n",
    "                \n",
    "    # 停止计时器并打印训练用时\n",
    "    end_time = perf_counter()\n",
    "    print(\"MBQC 运行时间为： \", end_time - start_time)\n",
    "\n",
    "    # 打印优化后的训练参数 gamma, beta\n",
    "    print(\"最优参数 gamma 为： \", mbqc_net.gamma.numpy())\n",
    "    print(\"最优参数 beta 为： \", mbqc_net.beta.numpy())\n",
    "\n",
    "    # 解码原问题的答案\n",
    "    solution_str = get_solution_string(vector_out, SHOTS)\n",
    "\n",
    "    # 将最优解带回目标函数，求出对应的最优值\n",
    "    relation = {symbols('x_' + str(j + 1)): int(solution_str[j]) for j in range(var_num)}\n",
    "    value = poly_symbol.evalf(subs=relation)\n",
    "    \n",
    "    # 返回最优解和最优值\n",
    "    opt = [solution_str, value]\n",
    "    return opt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "为了检验训练结果的正确性，我们在 `pubo` 中定义了 `brute_force_search` 函数，用遍历的方法在解空间里暴力搜索 PUBO 问题的解，方便与 MBQC 模型训练得到的解进行对比。\n",
    "\n",
    "```python\n",
    "from paddle_quantum.mbqc.QAOA.pubo import brute_force_search\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Main 函数\n",
    "\n",
    "主函数定义后，我们就可以输入参数运行啦！"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# 定义主函数\n",
    "def main():\n",
    "    \n",
    "    # 以上面的例子 x_1 + x_2 - x_3 + x_1*x_2 -x_1*x_2*x_3 + 0.5 为例，求解 PUBO 问题\n",
    "    var_num = 3  \n",
    "    poly_dict = {'x_1': 1, 'x_2': 1, 'x_3': -1, 'x_1,x_2': 1, 'x_1,x_2,x_3': -1, 'cons':0.5}\n",
    "    polynomial = [var_num, poly_dict]\n",
    "    \n",
    "    # 打印输入的多项式\n",
    "    print(\"输入的多项式为：\", polynomial)\n",
    "    \n",
    "    # 我们也可以随机生成一个目标多项式\n",
    "    # polynomial = random_poly(var_num)\n",
    "\n",
    "    # 检查输入多项式的合法性\n",
    "    is_poly_valid(polynomial)\n",
    "\n",
    "    # 进行训练并返回最优结果\n",
    "    mbqc_result = mbqc_pubo(\n",
    "        OBJ_POLY=polynomial,  # 目标多项式函数\n",
    "        DEPTH=6,  # QAOA 算法电路深度\n",
    "        SEED=1024,  # 随机种子\n",
    "        LR=0.1,  # 学习率\n",
    "        ITR=120,  # 训练迭代次数\n",
    "        EPOCH=1  # 迭代周期\n",
    "    )\n",
    "\n",
    "    # 打印 MBQC 模型训练出来的最优结果\n",
    "    print(\"MBQC 模型得到的最优解为：\", mbqc_result[0])\n",
    "    print(\"MBQC 模型得到的最优值为：\", mbqc_result[1])\n",
    "    \n",
    "    # 通过暴力搜索计算 PUBO 问题的最优结果并打印\n",
    "    brute_result = brute_force_search(polynomial)\n",
    "    print(\"暴力搜索得到的最优解为：\", brute_result[0])\n",
    "    print(\"暴力搜索得到的最优值为：\", brute_result[1])\n",
    "    \n",
    "    # 将两种结果进行比较\n",
    "    print(\"MBQC 模型和暴力搜索的差值为：\", mbqc_result[1] - brute_result[1])\n",
    "\n",
    "\n",
    "if __name__ == '__main__':\n",
    "    main()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "## 示例：最大割问题\n",
    "\n",
    "### 图与割\n",
    "\n",
    "最大割问题（MaxCut Problem）是图论中常见的一个组合优化问题，在统计物理学和电路设计中都有重要应用。\n",
    "\n",
    "在图论中，一个图由 $G = (V, E)$ 表示，其中 $V$ 为图的点集合，$E$ 为边集合。例如一个正方形就可以由 $G = (V,E)$ 其中 $V = [1,2,3,4]$， $E = [(1,2),(2,3),(3,4),(1,4)]$ 来表示。\n",
    "\n",
    "代码上，我们可以用 `maxcut` 中的 `plot_graph` 画出给定的图。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "from paddle_quantum.mbqc.QAOA.maxcut import plot_graph\n",
    "V = [1,2,3,4]\n",
    "E = [(1,2),(2,3),(3,4),(1,4)]\n",
    "G = [V, E] \n",
    "plot_graph(G,\"A square\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "一个图的割是指将该图的顶点集 $V$ 分割成两个互不相交的子集合 $S_0$ 和 $S_1$ 的一种划分。如果图中边的两个顶点被划分到不同集合中，那么记录得一分，总得分则定义为这个割的大小。最大割问题就是要找到一个割使得该割的大小最大。\n",
    "\n",
    "比如对应如上所述正方形图 $G$， 其中一个最大割就是将点 $1$ 和 $3$ 分为同一个集合 $S_0$ 中，点 $2$ 和 $4$ 分到另一个集合 $S_1$ 中。  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 转化为 PUBO 问题\n",
    "\n",
    "最大割问题其实可以转化为一个 PUBO 问题。假设待分割的图 $G = (V, E)$ 有 $n=|V|$ 个顶点和 $m =|E|$ 条边，那么我们可以将最大割问题等价为拥有 $n$ 个自变量的组合优化问题。每个变量 $x_v$ 对应图 $G$ 中的一个顶点 $v \\in V$，其取值 $x_v \\in \\{0,1\\}$，分别对应于顶点属于集合 $S_0$ 或 $S_1$，因此 $x = \\{x_1,\\cdots,x_n\\}$ 的每一种取值方式都对应一个割。由于最大割问题中有效计数的边为那些顶点属于不同集合的边，即顶点 $u$ 和 $v$ 对应的变量取值不同 $x_u \\neq x_v$ 才会被记一分。所以对于给定的割 $x$，其对应的割的大小为\n",
    "\n",
    "$$\n",
    "C(x) = \\sum_{(u,v) \\in E} (x_u \\oplus x_v),\\tag{3}\n",
    "$$\n",
    "\n",
    "其中 $\\oplus$ 为异或运算。那么最大割问题等价于求解 $ \\underset{x}{\\max} \\ C(x)$。进一步，上述函数 $C(x)$ 可以写成多项式形式\n",
    "\n",
    "$$\n",
    "C(x) = \\sum_{(u, v) \\in E} (x_u + x_v - 2 x_u x_v),\\tag{4}\n",
    "$$\n",
    "\n",
    "从而最大割问题等价于一个 $n$ 元 $2$ 次多项式的 PUBO 问题。我们希望找到最优解 $x^{*}$ 使得目标多项式最大，\n",
    "\n",
    "$$\n",
    "x^* = \\underset{x}{\\text{argmax}} \\ \\left( \\sum_{(u, v) \\in E} (x_u + x_v - 2 x_u x_v) \\right).\\tag{5}\n",
    "$$\n",
    "\n",
    "在 `maxcut` 当中，我们定义了 `graph_to_poly` 函数，可以输入待分割的图，返回等价的 PUBO 问题的目标多项式。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 引入 maxcut 模块\n",
    "from paddle_quantum.mbqc.QAOA.maxcut import graph_to_poly\n",
    "\n",
    "# 输入点和边的信息\n",
    "V = [1,2,3,4]\n",
    "E = [(1,2),(2,3),(3,4),(1,4)]\n",
    "# 构造待分割的图\n",
    "G = [V, E] \n",
    "\n",
    "# 将图转化为对应 PUBO 问题的目标多项式\n",
    "poly = graph_to_poly(G)\n",
    "print(\"与图等价的目标多项式为：\\n\", poly)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 求解最大割问题的代码实现\n",
    "\n",
    "得到等价的目标多项式之后，我们便可以仿照 PUBO 的流程，把最大割问题作为 PUBO 问题来求解。具体代码如下："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 引入符号运算模块\n",
    "from sympy import symbols\n",
    "# 引入 paddle 模块\n",
    "from paddle import seed, optimizer\n",
    "# 引入 qaoa 模块\n",
    "from paddle_quantum.mbqc.QAOA.qaoa import MBQC_QAOA_Net, get_solution_string\n",
    "# 引入 maxcut 模块\n",
    "from paddle_quantum.mbqc.QAOA. maxcut import plot_graph, graph_to_poly, plot_solution"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "利用 MB-QAOA 算法，我们可以定义 MaxCut 主函数，将图输入到 MB-QAOA 算法中，最终返回求得的解并画出割的方式。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 定义 MaxCut 主函数\n",
    "def mbqc_maxcut(GRAPH, DEPTH, SEED, LR, ITR, EPOCH, SHOTS=1024):\n",
    "    \n",
    "    # 打印待分割的图\n",
    "    plot_graph(graph=GRAPH, title=\"Graph to be cut\")\n",
    "\n",
    "    # 找到图对应的目标多项式\n",
    "    polynomial = graph_to_poly(GRAPH)\n",
    "    print(\"与图等价的目标多项式为：\", polynomial[1])\n",
    "\n",
    "    # 开始算法的记时\n",
    "    start_time = perf_counter()\n",
    "\n",
    "    # 实例化一个 MB-QAOA 训练网络\n",
    "    seed(SEED)\n",
    "    mbqc_net = MBQC_QAOA_Net(DEPTH)\n",
    "    # 选择 Adams 优化器\n",
    "    opt = optimizer.Adam(learning_rate=LR, parameters=mbqc_net.parameters())\n",
    "\n",
    "    # 开始训练\n",
    "    for epoch in range(EPOCH):\n",
    "        # 每次迭代后更新训练参数\n",
    "        for itr in range(1, ITR + 1):\n",
    "            # 训练参数并返回损失函数值和当前输出量子态\n",
    "            loss, state_out = mbqc_net(poly=polynomial)\n",
    "            # 根据损失函数值优化训练参数\n",
    "            loss.backward()\n",
    "            opt.minimize(loss)\n",
    "            opt.clear_grad()\n",
    "            if itr % 10 == 0:\n",
    "                print(\"iter:\", itr, \"  loss_MBQC:\", \"%.4f\" % loss.numpy())\n",
    "\n",
    "    # 停止算法的记时并打印用时\n",
    "    end_time = perf_counter()\n",
    "    print(\"MBQC 运行时间为：\", end_time - start_time)\n",
    "    \n",
    "    # 打印训练得到的最优参数\n",
    "    print(\"最优参数 gamma 为： \", mbqc_net.gamma.numpy())\n",
    "    print(\"最优参数 beta 为：\", mbqc_net.beta.numpy())\n",
    "\n",
    "    # 从量子态中解码问题的答案\n",
    "    mbqc_solution = get_solution_string(state_out, SHOTS)\n",
    "    # 绘制割的方式\n",
    "    plot_solution(GRAPH, mbqc_solution)\n",
    "\n",
    "    # 将最优解带入目标函数，求得对应最优值，并返回\n",
    "    var_num, poly_symbol = polynomial\n",
    "    relation = {symbols('x_' + str(j + 1)): int(mbqc_solution[j]) for j in range(var_num)}\n",
    "    \n",
    "    mbqc_value = int(poly_symbol.evalf(subs=relation))\n",
    "    mbqc_opt = [mbqc_solution, mbqc_value]\n",
    "\n",
    "    # 返回最优解和对应的最优值\n",
    "    return mbqc_opt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Main 函数\n",
    "\n",
    "主函数定义后，我们就可以输入参数运行啦！"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "def main():\n",
    "    # 以正方形的 MaxCut 问题为例，输入节点和边信息\n",
    "    V = [1, 2, 3, 4]\n",
    "    E = [(1, 2), (2, 3), (3, 4), (4, 1)]\n",
    "    G = [V, E]\n",
    "    \n",
    "    # 进行训练并返回最优结果\n",
    "    mbqc_result = mbqc_maxcut(\n",
    "        GRAPH=G,  #  待分割的图\n",
    "        DEPTH=6,  # QAOA 算法电路深度\n",
    "        SEED=1024,  # 随机种子\n",
    "        LR=0.1,  # 学习率\n",
    "        ITR=120,  # 训练迭代次数\n",
    "        EPOCH=1,  # 迭代周期\n",
    "        SHOTS=1024  # 解码答案时的采样数\n",
    "    )\n",
    "\n",
    "    # 打印训练的最优结果\n",
    "    print(\"最优解为：\", mbqc_result[0])\n",
    "    print(\"最优值为： \", mbqc_result[1])\n",
    "\n",
    "if __name__ == '__main__':\n",
    "    main()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "至此，我们完成了本教程的两个完整示例。MB-QAOA 算法实现预示着 MBQC 模型在量子机器学习领域强大的潜力。很显然，MBQC 模型所能够处理的算法不仅有 QAOA，我们非常期待 MBQC 这种非同寻常的计算方式能在某些特殊情境下展现出无与伦比的优势！"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "---\n",
    "\n",
    "## 参考文献\n",
    "\n",
    "[1] Farhi, Edward, et al. \"A quantum approximate optimization algorithm.\" [arXiv preprint arXiv:1411.4028 (2014).](https://arxiv.org/abs/1411.4028)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.10"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": false,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
